﻿using System;
using System.Text;
using System.Drawing;
using System.Buffers;
using System.Collections;
using System.Collections.Generic;
using System.Runtime.InteropServices;

public static partial class NativeAOT
{
    [UnmanagedCallersOnly(EntryPoint = "rhis")]
    public static unsafe void rhis(int n, IntPtr x_ptr, int m, double x0, double h, IntPtr g_ptr, IntPtr q_ptr, IntPtr dt_ptr, int k)
    {
        double* x = (double*)x_ptr.ToPointer();
        int* g = (int*)g_ptr.ToPointer();
        int* q = (int*)q_ptr.ToPointer();
        double* dt = (double*)dt_ptr.ToPointer();

        rhis(n, x, m, x0, h, g, q, dt, k);
    }

    /// <summary>
    /// 随机样本分析.
    ///  void rhis(int n, double x[], int m, double x0, double h, int g[], int q[], double dt[3], int k)
    /// </summary>
    /// <param name="n">随机样本点数。</param>
    /// <param name="x">x[n]存放随机变量的n个样本点值。</param>
    /// <param name="m">直方图中区间总数。</param>
    /// <param name="x0">直方图中随机变量的起始值。</param>
    /// <param name="h">直方图中随机变量等区间长度值。</param>
    /// <param name="g">g[m]返回m个区间的按高斯分布所应有的近似理论样本点数。</param>
    /// <param name="q">q[m]返回落在m个区间中每一个区间上的随机样本实际点数。</param>
    /// <param name="dt">dt[3]dt[0]返回随机样本的算术平均值，dt[1]返回随机样本的方差，dt[2]返回随机样本的标准差。</param>
    /// <param name="k">k=0表示不需要输出直方图；k=1表示需要输出直方图。</param>
    public static unsafe void rhis(int n, double* x, int m, double x0, double h, int* g, int* q, double* dt, int k = 0)
    {
        int i, j;//, kk, z;
        double s;

        dt[0] = 0.0;
        // 随机样本的算术平均值
        for (i = 0; i <= n - 1; i++)
        {
            dt[0] = dt[0] + x[i] / n;
        }
        dt[1] = 0.0;
        for (i = 0; i <= n - 1; i++)
        {
            dt[1] = dt[1] + (x[i] - dt[0]) * (x[i] - dt[0]);
        }
        // 随机样本的方差
        dt[1] = dt[1] / n;
        // 随机样本的标准差
        dt[2] = Math.Sqrt(dt[1]);

        // 按高斯分布所应有的近似理论样本点数
        for (i = 0; i <= m - 1; i++)
        {
            q[i] = 0;
            s = x0 + (i + 0.5) * h - dt[0];
            s = Math.Exp(-s * s / (2.0 * dt[1]));
            g[i] = (int)(n * s * h / (dt[2] * 2.5066));
        }
        s = x0 + m * h;
        // 落在每一个区间上的随机样本实际点数
        for (i = 0; i <= n - 1; i++)
        {
            if ((x[i] - x0) >= 0.0)
            {
                if ((s - x[i]) >= 0.0)
                {
                    j = (int)((x[i] - x0) / h);
                    q[j] = q[j] + 1;
                }
            }
        }
        //不需要输出直方图
        if (k == 0) return;

#if __UNUSED__
        // cout <<"n = " <<n <<endl;
        // cout <<"随机变量起始值 x0 = " <<x0 <<endl;
        // cout <<"随机变量区间长度 h = " <<h <<endl;
        // cout <<"直方图中区间总数 m = " <<m <<endl;
        // cout <<"样本算术平均值 = " <<dt[0] <<endl;
        // cout <<"样本的方差 = " <<dt[1] <<endl;
        // cout <<"样本的标准差 = " <<dt[2] <<endl;
        char[] a = new char[50];

        kk = 1; z = 0;
        for (i = 0; i <= m - 1; i++)
        {
            if (q[i] > z) z = q[i];
        }
        // kk为比例系数
        while (z > 50)
        {
            z = z / 2;
            kk = 2 * kk;
        }
        // cout <<"区间中点  实际点数    直方图" <<endl;
        for (i = 0; i <= m - 1; i++)
        {
            //区间中点值
            s = x0 + (i + 0.5) * h;
            for (j = 0; j <= 49; j++) a[j] = ' ';
            j = q[i] / kk;
            //实际点位置符号
            for (z = 0; z <= j - 1; z++) a[z] = 'X';
            j = g[i] / kk;
            //理论点数位置符号
            if ((j > 0) && (j < 50)) a[j] = '*';
            // cout <<setw(8) <<s <<setw(10) <<q[i] <<"   ";//
            for (j = 0; j <= 49; j++)
            {
                //cout << a[j];
            }
            // cout <<endl;
        }
        // cout <<"比例  1 : " <<kk <<endl;      
        return;
#endif
    }

    /*
    // 随机样本分析例1
      int main()
      { 
          int n,m,k,g[10],q[10];
          double dt[3],x0,h;
          double x[100]={
                  193.199,195.673,195.757,196.051,196.092,196.596,
                  196.579,196.763,196.847,197.267,197.392,197.477,
                  198.189,193.850,198.944,199.070,199.111,199.153,
                  199.237,199.698,199.572,199.614,199.824,199.908,
                  200.188,200.160,200.243,200.285,200.453,200.704,
                  200.746,200.830,200.872,200.914,200.956,200.998,
                  200.998,201.123,201.208,201.333,201.375,201.543,
                  201.543,201.584,201.711,201.878,201.919,202.004,
                  202.004,202.088,202.172,202.172,202.297,202.339,
                  202.381,202.507,202.591,202.716,202.633,202.884,
                  203.051,203.052,203.094,203.094,203.177,203.178,
                  203.219,203.764,203.765,203.848,203.890,203.974,
                  204.184,204.267,204.352,204.352,204.729,205.106,
                  205.148,205.231,205.357,205.400,205.483,206.070,
                  206.112,206.154,206.155,206.615,206.657,206.993,
                  207.243,207.621,208.124,208.375,208.502,208.628,
                  208.670,208.711,210.012,211.394};
          n=100; m=10; x0=192.0; h=2.0; k=1;
          rhis(n,x,m,x0,h,g,q,dt,k);
          return 0;

      }
    */
    /*
    // 随机样本分析例2
      int main()
      { 
          int n,m,k,j,g[10],q[10];
          double dt[3],x0,h;
          double x[500];
          RND p;
          p=RND(1.0);
        // 产生500个均值为100方差为2.25的正态分布随机数
          for (j=0; j<500; j++)
              x[j] = p.rndg(100.0, 1.5);
          n=500; m=10; x0=91.0; h=2.0; k=1;
          rhis(n,x,m,x0,h,g,q,dt,k);
          return 0;

      }
    */
}

